#ifndef EQUATIONSOLVER
#define EQUATIONSOLVER

#include "Function.h"
#include <iostream>

class EquationSolver {
protected:
    const Function& F;
public:
    EquationSolver(const Function& F) : F(F) {}
    virtual double solve() = 0;
    ~EquationSolver() {}
};

class Bisection_Method : public EquationSolver {
private:
    double a, b;
    double eps, delta;
    int Maxiter;
public:
    Bisection_Method(const Function& F, double a, double b,
        double eps = 1e-7, double delta = 1e-6, int Maxiter = 50) :
        EquationSolver(F), a(a), b(b), eps(eps), delta(delta), Maxiter(Maxiter) {}
    ~Bisection_Method() {}
    virtual double solve() {
        if (F(a) == 0) return a;
        if (F(b) == 0) return b;
        if (F(a) * F(b) > 0) 
        {
            std::cerr << "the input a and b have the same sign" << std::endl;
            return NAN;
        }
        double h = b - a;
        double u = a + h / 2;;
        double a1 = a;
        int k=0;
        for(int i=0;i<=Maxiter;i++)
        {
            k = i;
            u = a1 + h / 2;
            h = h / 2;
            if (h <= eps || abs(F(u)) <=eps ) return u;
            if (F(u) * F(a1) >0) a1 = u;

        }
        if (k == Maxiter)std::cout << "Reache the upper bound of the number of cycles" << std::endl;
        return u; 
    }
};

class Newton_Method : public EquationSolver {
private:
    double x0;
    double eps;
   double  Maxiter;
public:
    Newton_Method(const Function& F, double x0,
        double eps = 1e-7, double Maxiter = 8) :
        EquationSolver(F), x0(x0), eps(eps) ,Maxiter(Maxiter) {}

    virtual double solve() 
    {
        int k=0 ;
        double u;
        double x=x0;

        for (int i = 1; i <= Maxiter; i++) 
        {
            u = F(x);
            if (abs(u) <= eps) return x;
            if (abs(F.derivative(x)) <= eps) 
            {
                std::cerr << "derivative is too small ,no solution" << std::endl;
                return NAN;
            }
            x = x - u / F.derivative(x);
            k = i;
        }
        if (k == Maxiter)std::cout << "Reache the upper bound of the number of cycles" << std::endl;
        return x; 
    }
};

class Secant_Method :public EquationSolver
{
public:
private:
    double x0, x1;
    double eps, delta;
    int Maxiter;
public:
    Secant_Method(const Function& F, double x0, double x1,
        double eps = 1e-7, double delta = 1e-6, int Maxiter = 15) :
        EquationSolver(F), x0(x0), x1(x1), eps(eps), delta(delta), Maxiter(Maxiter) {}
    virtual double solve() 
    {
        double y0 = x0;
        double y1 = x1;
        double u = F(y0);
        double v = F(y1);
        int k = 1;
        if (abs(y1 - y0) <= eps || abs(v) <= delta) return y1;
        for(int i=2;i<=Maxiter;i++)
        {
            k = i;
            if (fabs(v - u) < eps) 
            {
                std::cerr << "Error: Division by zero in Secant method." << std::endl;
                return NAN;
            }
            double s = (y1 - y0) / (v - u);
            u = v;
            y0 = y1;
            y1 = y0 - v * s;
            v = F(y1);
            if (abs(y1 - y0) <= eps || abs(v) <= delta) break;
    
        }
        if (k == Maxiter)std::cout << "reach the upper bound of number of loops" << std::endl;
        return y1;
    }
};

#endif